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We propose a new method for detecting paired states in either bosonic or fermionic systems 
using interference experiments with independent or weakly coupled low dimensional systems. We 
demonstrate that our method can be used to detect both the FFLO and the d-wave paired states 
of fermions, as well as quasicondensates of singlet pairs for polar F=l atoms in two dimensional 
systems. We discuss how this method can be used to perform phase-sensitive determination of the 
symmetry of the pairing amplitude. 



I. INTRODUCTION 



Interference experiments are the primary tool of de- 
tecting and characterizing cold atom systems^ 2 -. While 
original experiments focused on demonstrating macro- 
scopic coherence of large BEC's jjj], subsequent work 
used interference experiments to explore more interesting 
phases and phenomena. For example, interference in the 
time of flight (TOF) experiments was used for observa- 
tion of the superfluid to Mott insulator transition in op- 
tical lattices 4 , analysis of fluctuations in low dimensional 
systems^, and studies of phase diffusion and decoher- 
ence in dynamically split condensates 3 -^. Interference 
can also give rise to interesting patterns in second order 
coherence^. This approach was used to demonstrate that 
Hanburry Brown Twiss experiments with both bosons 
and fermion o 10 ' 11 ' 12 ' 13 ' 1 ' 15 ' 16 ' and to observe pairing of 
fermions 17 . A series of recent theoretical and experimen- 
tal papers explored the idea that one can use interfer- 
ence between two or more low-dimensional systems to 
probe their non-trivial correlation functions 5 -! 6 -! 1 ^. Par- 
tially using this ideas Hadzibabic et. al. were able to 
detect Berczinskii-Kostcrlits-Thouless transition in two- 
dimensional bosonic systems which is associated with 
vortex proliferation 5 -. One of the novel feature of this 
approach was the idea to use not only the average con- 
trast but the full distribution functions 6 - ' 19 ' 20 . Distri- 
bution functions are determined by high order correla- 
tion functions and contain a wealth of information about 
underlying systems. Distribution functions of interfer- 
ence fringe amplitudes were recently analyzed for one- 
dimensional quasi-condensates and provided direct probe 
of long wavelength phase fluctuations in the system of ei- 
ther quantum or thermal origin [6| . However there is an- 
other source of fluctuations of the fringe amplitude which 
is purely quantum in nature. Namely, this is shot noise 
coming from the discreteness of particles 6 ^. Shot noise 
is especially strong in systems with short range single 
particle correlations, in particular in fermionic systems. 
Thus interferometric probes in such systems are intrin- 
sically more difficult than in the systems with long or 
quasi long range order for which shot noise is less impor- 



tant than the low wavelength thermal and/or quantum 
fluctuation s 20 ' 21 . Our emphasis on low-dimensional sys- 
tems has two main reasons: they exhibit exotic phases 
more often and it is easy to perform interference experi- 
ments with them. In this paper we focus on fermionic and 
bosonic systems in low dimensions and continue to study 
the possibility of using interferometry to probe strongly- 
correlated many-body states. 

In interacting systems one is often interested in states 
which do not have coherence of individual particles but 
exhibit a coherence (or slowly decaying correlations) of 
particle pairs. For example, fermionic paired states are 
characterized by the pairing amplitude 

A(r) = J drj f{r,) Cl (r - V /2) c T (r + rj/2) (1) 

Here r is the center of mass position of Cooper pairs 
and f(rj) is the Cooper pair wave function^ 3 -. The or- 
dered state corresponds to the condensation of pairs of 
particles and should be analyzed using correlation func- 
tions of the form (AT(n)A(r2)). Correlation functions of 
this type which we will refer to as anomalous correlation 
functions also arise in the context of exotic states of in- 
teracting bosons such as condensates of pairs of bosons 2 ^ 
and polar condensates in two dimensional system s 25 ' 26 . 
In principle one can extract anomalous correlation func- 
tions analyzing higher order moments of the interference 
amplitude. However, as we will show below, this might 
be a very difficult task in practice because of effects of 
shot nois o 20 ' 21 and because such anomalous correlation 
functions can appear as small corrections on top of nor- 
mal correlation functions. 

In this paper we suggest an alternative method for 
identifying paired states and for measuring directly their 
anomalous correlation functions using interference exper- 
iments with two (or more) systems. This paper extends 
earlier work on the analysis of interference experiments 
with pairs of independent condensates of single compo- 
nent bosons 5 -! 6 - ' 18 ' 19 ' 20 . Our main purpose here is to show 
that one can probe fermionic superfluidity in low dimen- 
sional systems. In particular, we define a new observable, 
which we refer to in the text as anomalous interference 
amplitude, which should vanish when there is no pair- 
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ing between fermions and which is nonzero when there is 
paring in the system. We suggest two methods to detect 
this anomalous amplitude. The first approach relies on 
detecting interferometric signal in two disjoint parts of 
the system Riji and averaging appropriate observable 
over these disjoint regions. This way of detecting pairing 
correlations relies on the existence of the long range (or 
quasi long range) order in the pairing channel and thus 
requires phase ordering in the fermionic superfluids. Note 
that averaging over two disjoint regions is necessary to 
cancel the effects of an undefined relative phase of the su- 
pcrfluid order parameter in two independent layers. One 
can straightforwardly extend this idea and split the sys- 
tem to a larger number of disjoint regions improving the 
signal to noise ratio but other than that not affecting 
our analysis. In the second method we introduce a weak 
tunneling coupling between the systems to lock the rel- 
ative phase. We show that in sufficiently large systems 
there is always a broad range of parameters, where the 
coherence is established but the correlation functions are 
still not affected by the presence of this weak tunneling 
term. Because the phase locking transition does not re- 
quire long range order in each superfluid, this method is 
more sensitive to the formation of the local pairing ampli- 
tude. We further argue that in lattice fermionic systems 
one can measure the symmetry of the pairing gap and 
thus distinguish, for example, d— wave from s— wave su- 
perfluidity This can be achieved by aligning the probing 
laser beam along different axes of the lattice. 

The ideas presented in this paper can be further ex- 
tended to low-dimensional Bose systems. We show that 
in a similar setup one can measure anomalous correlators 
in bosonic superfluids. These correlators have an unusual 
property that they grow with the separation between the 
particles showing effective "anti-bunching" behavior for 
bosons. Usually anomalous correlations are not easy to 
detect, since they are not gauge invariant, i.e. they are 
sensitive to the global superfluid phase. The two setups 
considered here eliminate effects of this phase and make 
such measurements possible. 

Carusotto and Castin have previously suggested an ex- 
periment which relies on particle interference to detect 
paired states^. While there is some conceptual connec- 
tion between their work and our approach, our method 
has an advantage that it does not require Bragg out- 
coupling of atoms, splitting and mixing of atom beams, 
and using single atom detectors to measure coincidences. 
As we demonstrate below, interference of two ballistically 
expanding independent clouds does all of this work itself! 

The paper is organized as follows. In Sec. [Til we 
first analyze the basic structure of anomalous correlators 
and the interference amplitude between two independent 
fermionic superfluids. We then introduce the new ob- 
servable, the anomalous interference amplitude, which 
probes the pairing amplitude. In Sec. IIHI we show how 
this anomalous amplitude can be detected performing si- 
multaneous measurements in disjoint parts of the time of 
flight image. Using this scheme we discuss possible set- 



ups for observing the d-wave superfluid and the FFLO 
phases. We suggest how one can detect not only the am- 
plitude, but also a phase of the pairing function. We 
perform explicit quantitative analysis of the anomalous 
amplitude for two-dimensional superfluids with s— and 
d— wave pairing based on BCS-theory. Then in Sec. IIVI 
we discuss the second way of detecting anomalous inter- 
ference amplitude by introducing a weak interlayer tun- 
neling. We show that on the one hand its presence in- 
troduces corrections to the results of Sec. |ITJ which are 
not related to the superfluidity. On the other hand the 
presence of this tunneling establishes the interlayer phase 
coherence. We show that by decreasing the imaging area 
and increasing the system size one can always achieve 
the regime where the coherence between the superflu- 
ids is established and yet the effect of the tunneling on 
the correlation functions is negligible. In Sec. [V] we ex- 
tend our analysis to bosonic superfluids. In particular, 
we show that in the superfluids with quasi long range 
order the anomalous interference amplitude grows super- 
linearly with the imaging size A. In turn this implies that 
the corresponding interference contrast increases with A. 
This behavior is opposite to that of the normal interfer- 
ence amplitude, which always decreases with A. And 
finally in Sec. I VII we summarize our results. 

Throughout the paper we use BCS approximation to 
perform explicit calculations. This approximation is only 
reliable in the weak coupling regime; at strong coupling 
one has to do more elaborate calculations. However, we 
do not expect any qualitative difference between BCS 
and exact results. 



II. ANALYSIS OF THE INTERFERENCE 
AMPLITUDE: BASIC SET-UP 

We start our discussion from analyzing the interfer- 
ence amplitude of two fermionic condensates. Extension 
of our results to the case of a stack of several condensates 
is straightforward. For concreteness we will focus on the 
case of two dimensions. First we analyze the usual inter- 
ference amplitude, which is determined by normal cor- 
relation functions and show that it is not a reliable de- 
tection tool of superfluidity. Then we describe how one 
can use the same interference experiments but analyze 
the results differently to extract anomalous correlation 
functions. 



A. Normal correlation functions. 

Consider two independent systems (layers) and assume 
that each system contains two species of atoms, which we 
label by a spin index a. Let c\ a (r) be the creation op- 
erators for atoms with spin a in layer i = 1,2 and the 
in plane coordinate r. After the expansion we find in- 
terference fringes in the z direction, so that the density 
Pint,o-(z,r) — Ccr(r) cos(Qz + 0<r(r)) where Q = md/Ht 
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(this assumes sufficiently long expansion time, see e.g. 
Ref. [H[). Because the phase </> CT (r) is a random vari- 
able for independent systems the average density does 
not show any interference fringes. Thus to filter out this 
oscillating component we have to consider Fourier trans- 
form of the density-density correlation function. Indeed 
one can chose the following operator, which corresponds 
to the square of the interference amplitude 2 ^: 



dz\dz2dv\dxi =F / p(z,r,t)dzdr, 



(2) 



where the upper (lower) sign corresponds to bosons 
(fermions). Here p{z, r, t) is the atomic density at posi- 
tion z, r at time t after the expansion. The z coordinate 
is orthogonal to the atomic systems, while r describes po- 
sitions of the atoms within each individual system. For 
two dimensional systems, integration over one of the di- 
rections is done automatically by the laser beam, whereas 
integration in the other direction is done manually^. We 
assume that the transverse confinement is tight and when 
the atoms are released, they expand strongly in the trans- 
verse z direction, while their in plane expansion can be 
neglected. This assumption is well justified if the trans- 
verse confining energy is large compared to any other 
energy scales in the problem. 

Before explaining where the expression ^ came from 
let us investigate it a little further. Assuming that the 
long time of flight allows us to use the far field expres- 
sion a 20 ' 21 i 22 we find 

\ A \ 2 =J2 j / dr l dr 2cL( I 'l)4, CT '( r 2)c2, CT '(l-l)ci^(r2). 

a, a' 

(3) 

Note that both for bosons and the fermions the expres- 
sion above can be obtained from the complex interference 
amplitude defined as 



A = Y,J dvA ^ v ) = j dvc\ a (v)c 2a {v). 



(4) 



On can think about A as of the Fourier transform of the 
density of the expanded cloud in the z direction (for a 
given r) 61 . Then the expression ([3|) can be obtained as 
the normal ordered product of A' A: 



\A\ 2 =±:A^A:, 



(5) 



where the "+" sign corresponds to bosons and the "— " 
sign does to fermions. 

In the case of independent systems there is no coher- 
ence between atoms hence expectation value of A is zero. 
This does not mean the absence of interference fringes in 
individual shots but only tells us about the random phase 
of interference fringes. Indeed the quantity \A\ 2 is insen- 
sitive to this phase. It directly measures the (square of 
the) amplitude of the interference and it does not average 
to zero even for independent systems. 



Let us make a few comments on where the expressions 
above come from. In the Eq. @ we are taking Fourier 
transform of the product of the densities of atoms af- 
ter expansion. This Fourier transform picks the compo- 
nent in this product oscillating with the wavevector Q 
and thus corresponding to the interference between the 
two systems. Note that the operator \A\ 2 in Eq. @ is 
very similar to the one originally introduced for bosons^ 
except for the negative sign appearing for fermions and 
except for the additional second term. The negative sign 
takes care of the fermionic statistics, or equivalently of 
the additional it phase shift in the interference part of 
the density-density correlation functions 2 ^. The second 
term in Eq. @ removes the trivial contribution to the 
Fourier transform coming from shot noise which is not 
related to the interference. This term is usually unim- 
portant for bosonic systems. Note that Eq. ([2]) can be 
rewritten using the normal ordered product of densities: 



\A\ 2 = J : p(z 1 ,r 1 ,t)p(z 2 ,r 2 ,t) : e K 



<{z 1 -z 2 ) 



(0) 



Substituting the far field expansion of the bosonic oper- 
into Eq. ^ we easily recover Eq. We 
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emphasize that it is important to first take the normal 
order in the product of densities p(z\, T\, t)p(z2, r 2 , t) and 
only after use the far field expansion for the density oper- 
ators. Using the opposite order will give spurious contri- 
butions. In bosonic systems with large number of atoms 
in the same state the creation and annihilation opera- 
tors can be approximately treated as commuting classi- 
cal fields and thus no ambiguity with ordering appears 
and the shot noise contribution is smal l 20 ' 21 . However, 
for fermionic systems, where the shot noise is usually im- 
portant one has to be careful in evaluating integrals like 
those appearing in Eq. d2|). 

For bosons statistical and scaling properties of \A\ 2 
contain important information about superfluidit y 18 ' 19 , 
which can be straightforwardly detected in experi- 
ments^. At the same time for fermions information 
about superfluidity is encoded in the Cooper pair cor- 
relation functions. Although the operator \A\ 2 certainly 
contains the information about superfluidity (see Ap- 
pendix [X| and can be in principle used to determine the 
pairing, it does not provide a "smoking gun" for detecting 
superfluidity. Indeed pairing only quantitatively affects 
the magnitude of the interference amplitude \A\. This 
magnitude can be affected also by various other reasons. 
Thus it is important to find another observable which 
vanishes unless fermions are paired. We are going to in- 
troduce such an observable in the next section. 



B. Anomalous correlation functions 

An observable, which directly probes the pairing wave 
function can be constructed from Eq. ^ with a slight 



modification: 

A Q = J p{zi^ ll t)p{z2 1 Y2,t)e lQ{zi+z ^dz 1 dz 2 dr l dv2. 

(7) 

Note the difference between Eqs. ^ and |J7J|. The former 
corresponds to taking the product of the Fourier trans- 
form of the density p{z) and its complex conjugate. The 
latter corresponds to taking the square of the Fourier 
amplitude without taking complex conjugation. For the 
long expansion time this expression reduces to 

A " = E / 7 dr^cljrM^i^) 

a, a' 

x c 2 , CT '(r 2 )c2, CT (ri). (8) 

The quantity A 2 =: A 2 : looks like exactly what we 
need. Indeed for independent layers it depends only on 
the product of the pairing amplitudes in the two layers: 

A 2 = (A 2 ) = J J dr 1 dr 2 F*(r 1) v 2 )F 2 {v 1 ,r 2 ), (9) 

where F*(ri,r 2 ) = (cL(ri)c|, (r 2 )). However, there is 
one subtlety. Unlike the normal amplitude squared \A\ 2 , 
which is always a positive real number, the anomalous 
amplitude squared A 2 is complex. Moreover for indepen- 
dent condensates A 2 is equal to zero because the phases 
of F% and F 2 are not correlated. To avoid this phase 
uncertainty one can try to look into |A 2 | 2 , which will 
involve second order correlation functions in each layer. 
However, it is easy to see that \A 2 \ 2 will be dominated 
by shot noise and normal (not anomalous) correlation 
functions. Thus there will be no advantage compared to 
analyzing A 2 . 

The main purpose of this paper is to show that one 
can overcome the effect of the uncertain relative phase 
and measure the anomalous amplitude A and thus detect 
superfluidity in fermionic systems. We note that funda- 
mental reason why extra efforts are needed to measure 
A 2 compared to A 2 is because the former (for indepen- 
dent systems) is not a gauge invariant quantity. This 
difficulty will similarly arise if one tries to measure not 
gauge invariant quantities in other setups. The ideas of 
this work can be extended to those situations as well. In 
later sections we will discuss some other examples of this 
kind. 

Here we suggest two different setups to fix the problem 
with the unknown phase in Eq. @. In the first setup, 
which we refer to as Scheme I (see Fig. (JTJ), we get rid of 
the random phase by making a special choice of the spa- 
tially separated integration domains. Namely instead of 
integrating A 2 over the entire region one splits the imag- 
ing area spanned by ri and r 2 to two spatially separated 
domains Ri and Rn. In each experimental run one inde- 
pendently determines A 2 in the two domains then takes 
the absolute value of their square and averages over many 
experimental runs. As we will show in detail this setup 
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FIG. 1: The experimental set-up corresponding to the Scheme 
I of our approach. The integration in X-direction is performed 
automatically by the imagine beam, whereas the integration 
in the Y-direction is done by " hands" . The interference signal 
is collected from two disjoint regions Ri and Rh. 

relies on the fact that single particle correlation func- 
tions decay sufficiently fast with the distance, while the 
pair correlation functions decay slowly or do not decay 
at all. This setup has an obvious advantage compared to 
measuring \A 2 \ 2 because single particle normal correla- 
tion functions decay fast with the distance. As a result 
the quantity A 2 (Oi)j4* 2 (f2 2 ) is dominated by anomalous 
correlation functions: 

( c \,i c \j c ^ IlC l' n } ~ ( c T.-f c l,-f^ c b- rjC T.")' 

where subscripts I and II indicate that spatially these 
operators are located in regions Ri and Rij. Note that 
if the domains are not spatially separated or single par- 
ticle correlations functions do not decay fast there are 
additional (unwanted) cross correlations in the equation 
above like (cj 7 cj jj)(c{ jC^u). 

In the second setup, which we refer as Scheme II, (see 
Fig. @) we introduce a weak tunneling t± between the 
two layers. This tunneling locks the phases of pairing 
amplitudes in the two layers and makes the expectation 
value of A 2 real and positive. Besides the phase lock- 
ing effect, the tunneling induces the mixing between the 
fcrmions in the two layers and results to the nonzero con- 
tribution to A 2 in Eq. ([8]) even in the absence of pair- 
ing. Below we will show that at small temperatures 
it is always possible to choose the tunneling such that 
phase locking transition already occurred but the corre- 
lation functions are not yet significantly affected so that 
Eq. (112[) still holds. The two setups are complimentary 




FIG. 2: The experimental set-up corresponding to the Scheme 
II of our approach. The non-zero tunneling between two lay- 
ers locks the phases of superfluids. The interference signal is 
collected now from the whole areas of superfluids. 

to each other and can be used depending on the situa- 
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tion. The Scheme I essentially relies on the existence of 
the long (or quasi-long) range order in the pairing am- 
plitude. As we will show this scheme can be adopted 
to measuring not only the existence of superfluidity but 
also to the symmetry of the order parameter and even 
its phase. While the Scheme II is more sensitive to the 
local pairing between fermions and less to the existence 
of long range order in the superfluid phase. Such setup 
can be used, for example, to measure pseudo gap phe- 
nomena. Scheme II can also be used to determine the 
local symmetry of the pairing wave function but not its 
phase. 



III. SCHEME I: BASIC SET-UP AND VARIOUS 
EXAMPLES 

Keeping the analogy with analysis of normal correla- 
tions we emphasize the integration region R in the defi- 
nition of the operator A 2 and denote it as A 2 (R) in what 
follows 

A 2 {R I ) = A ] {R I )A i (R I ) = 
= / dri / dr2 4 T (ri)cl|(r 2 )c2i(r2)c2 T (ri)(10) 

JR! JR., 

Since for independent condensates the correlation func- 
tion above factorizes into a product of anomalous cor- 
relation functions in each system and from Eq. (|10[) we 
arrive to Eq. ([9]). As we argued earlier, because A 2 (Rj) 
is a complex number with a phase which is random from 
shot to shot, taking expectation value of equation © 
gives us zero. To get rid of this random phase we com- 
pare interference patterns from a pair of regions, Ri and 
Rji. More precisely we take [A 2 (Rj)Y A 2 (Rjj) , so that 
the random relative phase between pairing functions F\ 
and F 2 drops out. Experimentally this procedure cor- 
responds to taking a square of the Fourier transform of 
the density along x-direction integrated over the region 
I and multiplying it by a complex conjugate of a simi- 
lar quantity integrated over the region two. The result 
of this manipulation is then averaged over many experi- 
mental runs. Assuming that the system has a true long 
range order, taking regions Ri and Rh to be separated 
by a distance which is appreciably larger than the size of 
the Cooper pairs, and taking Rj and Rh to be identical 
we find 

([A 2 (R I )^A 2 (R U )) = 
/ dr 1 [ dr 2 F*(r 1 ,r 2 )F 2 (r 1 ,r 2 ) (11) 

JR! JR! 

We now consider several specific examples in which 
analysis of expectation values of the type (fTTj) can be 
used to identify interesting many-body states. 



A. Analysis of the anomalous interference 
amplitude within the BCS theory 

We now consider the integral in Eq. (|TT|) . To simplify 
calculations we also assume the translational invariance 
in both systems. In this case the pairing wave func- 
tion depends only in the difference between ri and r 2 : 
F Q (n,r 2 ) = F Q (n -r 2 ). Then 

A 2 =2A [ drFj(r)F 2 (r), (12) 

JR 

where the integration is again taken over the part of a 
condensate R] or Ru with the imaging area A. As we 
noted once the relative phase is taken care of and assum- 
ing the two condensates are identical we have F% = F 2 
and thus 

f d\c 

A 2 = 2Aj j^\F(k)\ 2 , (13) 

The integral above can be easily evaluated within the 
BCS model (we take zero temperature limit) 

A f dk 

A2 -2j(2^M- ™ 

where the pairing function Ak has to be specified for 
concrete type of pairing, = y/t£ + A£ and £k = 
k 2 /2m — /i is a single-particle dispersion. Thus if the 
pairing gap is isotropic and energy independent Ak = A 
then in 2D we find 

4*) = ^°arctan^°. (15) 

where the 2D constant density of states po is introduced. 
In the weak coupling limit Ao <C p, where p is a chemi- 
cal potential, we have A 2 ss 7VA§/4/i 2 , where N ps 2ppo 
is the total number of particles (the factor of two takes 
into account two different spin components). If the pair- 
ing is strong then BCS extrapolation gives A 2 ss Nir/4. 
We see that A 2 is a monotonically increasing function 
of the pairing gap and thus can serve as a direct probe 
of the fermionic superfluidity. Note that Eq. (fT"4"]) can 
be also analyzed in the case of d-wave pairing, where 
A(k) ~ Aocos(20k), where #k is the polar angle of the 
wave vector k. The result is (see Appendix C) 

However, since only the square of Ak enters Eq. (fl4|) 
and we are explicitly averaging over angles, the difference 
between s— and d— parings will be minor. In fact one can 
show that in the d— wave case Eq. (IT5|) gets multiplied by 
a smooth function of A a /p which changes between 1/2 
at Ao <C p and 2 /it in the opposite limit. 
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B. Phase sensitive detection of the d-wave pairing 

In the section above we discussed a possibility to 
detect anisotropy of the pairing amplitude using one- 
dimensional integration. In particular for the d-wave par- 
ing (dSF) the interference signal should vanish along the 
nodal directions. There are also other earlier suggestions 
for the detection of dSF, which rely rely on the detection 
of the Dirac like dispersion of quasiaprticle o 29 ' 30 ! 31 ' 32 . 
This however is not a unique signature of the d-wave 
pairing state. A Dirac cone of quasi-particles may also 
arise for an anisotropic s-wave pairing— or d-density wave 
states 3 -!. 

Here we would like to show how the Scheme I can be 
extended to do phase sensitive detection of dSF. In high 
temperature cuprate superconductors, the crucial exper- 
iments which identified the d-wave character of pairing 
were phase sensitive experiments by Van Harlingen et 
ali^ and Tsuei and Kirtley 3 ^ These experiments unam- 
biguously demonstrated the correct angular dependence 
of the pairing amplitude. Experimental set-up by Van 
Harlingen et al. used a combination of an s-wave and d- 
wave superconductors in a corner SQUID geometry. In- 
terference of s-wave Cooper pairs with different parts of 
d-wave Cooper pairs was used to establish the relative 
phase of of the Cooper pair wave function. 

What we discuss below is the cold atoms analogue of 
the Van Harlingen experiments. Hence we also need a 
source of s-wave Cooper pairs and a source of d-wave 
Cooper pairs. We imagine a pair of two dimensional 
fermionic systems, made of the same species of atoms, 
but having s-wave pairing in one layer and d-wave pair- 
ing in the other layer. This may be achieved, for example, 
using magnetic field dependence of the scattering length 
and applying a strong field gradient. Now we analyze in- 
terference patterns from two regions, i?j and Ru, which 
differ only by the 90° rotation. The quantity Aq(R]) 
is a complex number which has a random phase from 
one shot to another. Analogously Aq(Rji) is a complex 
number with a random phase. But the d-wave symmetry 
of the pairing requires that phases of these two complex 
amplitudes differ by precisely it. Hence one can look 
at ([Aq(Rj)]^ Aq(Rji)) and the d-wave symmetry dic- 
tates that this expectation value should be negative. On 
the other hand, when Rj and Rji have the same orien- 
tation, expectation value of (L4g(i?j)]^ Aq(Rji)) should 
be a positive number. It is important to emphasize that 
this statement is general and does not rely on the spe- 
cific microscopic model for d-wave pairing. We stress that 
that only one of the layers should have a d-wave symme- 
try otherwise Aq, being proportional to the product of 
anomalous correlation functions in two layers, does not 
change sign under rotations (see Eq. (|12jl ). While the pre- 
cise value of ([Aq(Ri)]^ Aq(Rh)} is not easy to calculate, 
especially if we are dealing with non-identical superfluids, 
that statement of the ir phase difference between Aq(Ri) 
and B2(Rii) relies only on the d-wave nature of pairing. 

The crucial feature of the method discussed in this sub- 



section is that it should provide a qualitative and model 
independent signatures of d-wave pairing. It does not 
rely on detailed analysis of the microscopic models but it 
uses only the fundamental symmetry of the d-wave order 
parameter. 



C. Probing of the anisotropy of pairing amplitude 

We now discuss another probe of d-wave pairing. Un- 
like the previous method, it can not be used to demon- 
strate the change of the sign of the gap function A(/c). 
However it can be used to observe anisotropy of the gap. 
We study the correlation function ([8]) but integrate it in 
a highly anisotropic way. In particular, one length, say 
along the probing beam should be macroscopic and the 
other W should be shorter than the coherence length. 
Then (the square of) the anomalous interference ampli- 
tude becomes 

a 2 {9) = 2AW J dz\F(z,8)\ 2 , (17) 

where 9 is the polar angle which defines direction of inte- 
gration. We introduced a new notation a 2 to avoid pos- 
sible confusion with A 2 analyzed earlier. Note that typ- 
ically d-wave symmetry of the order parameter requires 
the presence of the optical lattice. This lattice in turn 
breaks rotational symmetry in the superfluid and locks 
the phase of the pairing amplitude with the lattice's prin- 
cipal axes. Therefore there is no ambiguity in defining 9 
from one experimental run to another. One can expect 
that for s— wave pairing (fT7|) should give isotropic re- 
sult, while for the d— wave pairing the outcome will be 
highly anisotropic. While this approach does not provide 
a "smoking gun" signature of the change of sign in the 
Cooper pair wavefunction, this method is easier to do 
experimentally; if successful it should provide a strong 
indication of anisotropic pairing. 

It is straightforward to show that in the s— wave case 
the function F is isotropic and is given by 

F(z, 9) = Po A J (z/()K (z/£), (18) 

where 



hp = y / 2mp/h 2 and fc A = ^2mA/h 2 . This expression 
shows that the pairing wavefunction diverges logarithmi- 
cally at small z and decays exponentially with the char- 
acteristic correlation length £ at large z. The logarithmic 
divergence is the usual artifact of the BCS theory with 
point-like interactions. This divergence is cutoff at short 
distances. 

Using Eqs. (|18p we evaluate the integral in Eq. (fTOl) 
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and find 



a-2 



X 4^3 
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1 
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2' 


I) 



,(1,1,1),- 



c 2 



(19) 



where 4F3 is the generalized Hypergeometric function. 
At small and large ratio of A//j the expression above 
gives the following asymptotics: 



a 2 



«2 



Sir 
18.8 



,2 A 2 



AW 



In" [ * 



'nA 2 



/■'A 



A > /i 



A < w (20) 
(21) 



As before the high energy estimate of the asymptotics is 
the extrapolation of the BCS result to the strong cou- 
pling limit. For d~ wave pairing the correlation func- 
tion F(zg, 0) vanishes along the nodal direction and thus 
clq(0) should vanish as well. On the other hand, along 
the antinodal direction we can recover the asymptotics 
similar to s— wave case, 



a-2 



^AwEsg-iJ A <„. (22) 



For details of computations see Appendices B and C 



thought of as modifying the integral in equation J3|) 
A Q , a (R) = [ A a (r)e^ x dv = [ dve^ c\ a {v)c 2a {v) 

JR J R 

(23) 

where we implicitly assume that the direction of the 
vector Q coincides with the direction of integration x, 
Q = (<2,0). Defining now 

A 2 (R) = (A Q>a (R)A Q>a ,(R)) = 

= [ dx x ( dr 2 e lQ ^ 1+2 ' 2 )Fi t (r 1 ,r 2 )F 2 (r 1 ,r 2 ).(24) 

JR JR. 

In the FFLO phase A2 should have additional peaks at Q 
matching the finite momentum of the Cooper pair. The 
global unknown relative phase can be removed again ei- 
ther by multiplying the signal coming from two spatially 
separated imaging areas Ri and Ru or by introducing 
weak tunneling coupling between the layers as discussed 
in the next section. 

One may be concerned that in rotationally invariant 
systems the direction of the FFLO ordering wavevectors 
will not generically coincide with the x axis used for the 
observation. This issue should be avoided by using sys- 
tems that do not have a rotational symmetry in the xy 
plane. In fact, one of the most promising systems for 
observing the FFLO phase is an array of weakly coupled 
Id systems 40,41 . In this case the ordering wave vector 
should be in the direction of the tubes. 



D. FFLO phase 

One of the most intriguing suggestions for the paired 
states of fermions with attractive interactions is the idea 
of FFLO phase for systems with spin imbalance. This 
state corresponds to Cooper pairing at a finite momen- 
tum and has been a subject of extensive theoretical stud- 
ies during the last couple of year o 37 i 38 i 39 i 40 i 41 i 42 . Exper- 
imental situation remains unclear (for recent review see 
Ref. [Hj]). We now discuss how interference experiments 
can be adopted to look for the FFLO state. An ear- 
lier proposal for the detection of the FFLO phase can be 
found in Ref. [i3 |. 

The FFLO phase is characterized by the finite cen- 
ter of mass momentum of the Cooper pairs so that 
F(ri,r 2 ) = (c T (ri)q(r 2 )) oc e ±jQ ( ri+I " 2 )/ 2 . Therefore 
when one analyzes the anomalous interference amplitude 
one expects additional modulations which can be de- 
tected by taking an appropriate Fourier transform. This 
can be achieved by changing the integration procedure 
in Eqs. ([4]) and (fl2j) . We note that this integration is 
not equivalent for the two directions. In the direction 
of the x axis, integration is done automatically by the 
laser beam. In the other direction, i.e. along the y 
axis, it is performed "manually" by integrating interfer- 
ence fringes (see Fig. [2]). An alternative approach is to 
take a Fourier transform of the interference amplitude 
A a (y) = J dxA rT (x, y) along the y axis. This can also be 



IV. SCHEME II: ANOMALOUS CORRELATION 
FUNCTIONS IN PHASE LOCKED SYSTEMS 

Another way to overcome the effect of unknown rel- 
ative phase detecting anomalous interference amplitude 
is to introduce a weak tunneling between the two layers 
(see Fig. [SJ. As it was shown in Refs. [7f>. Hi such tun- 
neling leads to the phase-locking transition. At the same 
time if the tunneling is sufficiently weak then correlation 
functions do not appreciably change and Eq. (Tj"2"]) is still 
valid. Below we will show that there is indeed a wide 
range of parameters where the phases between pairing 
amplitudes in two layers are locked and Eq. (jT3J) gives 
the dominant contribution into the expression |8|). 

In the next section we analyze the effect of weak cou- 
pling more carefully. We will explicitly analyze only the 
case of two coupled s-wave super fluids. However, our 
results should be very general because precise nature of 
the symmetry of the pairing amplitude (s-wave, d-wave, 
FFLO, etc.) is not very important for the phase-locking 
phenomena. 

A. Role of the inter-layer coupling 

Two imaging areas used in the set-ups discussed previ- 
ously were needed to cancel the unknown relative phase 
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between the order parameters in two layers. The same 
effect, however, can be achieved by introducing a weak 
tunneling coupling between the layers. Then phases of or- 
der parameters should loc h 45 ' 46 and one does not have to 
combine signals from two different areas. Conversely dif- 
ferent imaging areas can be used as independent sources 
so that one can effectively average A2 over several inde- 
pendent imaging areas in a single experimental run. 

As before we will work in the BCS limit. The BCS 
Hamiltonian of two coupled condensates reads 

H = Y, (k) {euT* + A k r x ) Vi a (k) 

k. a 

-t± £>I(k)r,^(k) + v4(k)r 2 Vi(k)), (25) 



no need to go to r larger than the coherence length. The 
third term A 2 ^ comes from the Andreev process. It is 
always subdominant at small t± and we can safely ignore 
it. 

For the s— wave pairing the integrals above can be 
explicitly evaluated: 



A Po A A 

1 ~ arctan — , 

2 p 



A 



2,2 



At±p 



(31) 
(32) 



In the weak and strong pairing regimes we have the fol- 
lowing asymptotics: 



where we used Nambu notations: i\) a . \ = c a ,\, Vvi = 
c a a = 1,2 corresponds to two different layers, and 
t x and t z are the Pauli matrices. It is convenient to 
introduce symmetric and antisymmetric combinations: 
ip+ = + V'2)/v / 2 and = {ip x - ip 2 )/V2- The the 
Hamiltonian splits into the symmetric and antisymmetric 
parts H = Tt + + 7Y_: 

H± =^4(k)[(£ k ±i ± )r, + A k r,]V±(k). (26) 



Let us now find the effect of t± on A2 using Eq. ([8]). 
Reexpressing the operators c a and cj, through the Nambu 
spinors if) and ip' and using Wick's theorem and expand- 
ing to the leading order in the tunneling coupling t± it 
is straightforward to show that 



A 2 = ^2.1+^2,2+^2,3, 



where 



Ani fa — 



dk 



A k 



A 

2 ./ (2ir) d El 



A 



2.2 



!2,3 



AH\ 



2At 2 , 



dk 



(2ny El 
dk A 4 



(27) 

(28) 
(29) 
(30) 



A k - 



Here E£ = y/fa ± t^f 

Note that A%,x coincides with our earlier expression 
(fT4)l . The two other terms A 2 ,2 and A 2 ^ are proportional 
to We point that A 2 2 scales faster with the imaging 
area than the other terms. The reason is that even in the 
absence of superfluidity the tunneling forces fermions to 
occupy preferably the symmetric state. This is in turn 
equivalent to establishing a well defined relative phase 
between the two atomic systems. Since we are interested 
in detecting A 2 ^\ and not A 2 2 in fcrmionic systems, con- 
trary to bosonic, it is preferable to make the imaging area 
as small as possible, of the order of the Cooper pair size 
(or superconducting coherence length). Indeed the inte- 
gral in Eq. (fT2")) converges at long distances thus there is 



^^,AjKA*p*tl, A« M (33) 
Ai, 2 w ^A Po A, A\ « ~AH 2 ±P l A » Li. (34) 

Both at weak and strong pairing we find that 
A 2 . 2 



(35) 



where n is the atom density. The A 2) 2 contribution is 
an unwanted correction coming from the interlayer cou- 
pling, which have nothing to do with superfluidity. This 
contribution can be suppressed either by decreasing the 
tunneling amplitude t ± or by decreasing the imaging area 
A. Note that the tunneling amplitude can not be pushed 
down all the way to zero, because then one will loose 
phase coherence between the two layers. In the next sec- 
tion we will see that it is always possible to find the 
regime where ^2,2 is negligible and at the same time the 
two superfluids are locked. We also comment that one 
can distinguish two contributions A2.1 and ^2,2 by look- 
ing into the dependence of A2 on t±. 



B. Phase locking transition 

In this section we examine the effect of the tunneling 
coupling on establishing phase coherence between the two 
fcrmionic superfluids. From a prior work it is known 
that tow energy phase fluctuations in superfluids can be 
described by means of a conventional xy model with the 
effective Lagrangia n 47 ' 48 ! 49 : 



C = n6 



where 



- — (V6>) 2 - — Vn 2 - V - o; 

4m v ; I6mn y ' 2p 



n = n Q - Po (9+{V9 2 )/4m) 



(36) 



(37) 



plays the role of the superfluid density, 9 is the phase 
of the superfluid order parameter. We emphasize that 
if fluctuations of 9 are small then n is close to the total 
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density of atoms. It is straightforward to generalize the 
derivation of the Lagrangian given in Ref. [48[ to the 
bilayer system. A similar derivation but for a tunneling 
through a point can be also found in Ref. [47| and it is 
sketched in Appendix B for completeness. We use this 
formalism to evaluate the tunneling Hamiltonian from 
the effective action given in Eq. (|C9[) . 

In the case of s— and d-wave pairings the integral ap- 
pearing in Eq. (|C9|) can be easily evaluated yielding the 
following coupling term to the Hamiltonian of the system: 

t 2 

Hia = —±poT(A/vl) cos{9 1 - 6 2 ). (38) 

where the functions T(A//i) for the s— and d— wave pair- 
ings are given in Appendix C. In both cases T(A//x) is a 
smooth function, which interpolates between T(A//i) « 
2 at A <C p, and T(A//i) w 1 in the opposite limit. When 
the pairing is small A <§; // the effective Josephson cou- 
pling becomes simply J « t\po, he. independent on A. 
This is a bit surprising result since one would naively 
expect that J should vanish as A — * 0. And this is 
indeed the case in superfluids with inhomogeneous tun- 
neling where J ~ ^iPoA (see Ref. [47[ for details). In 
superfluids with uniform tj_ the Josephson coupling gets 
enhanced by the coherence factor l/(poA) ~ where 
£ is the BCS coherence length, and the dependence of J 
on A disappears. 

The minimal tunneling required to lock the two phases 
together between two superfluids can be estimated from 
equating the energy gap required to transfer one particle 
from one superfluid to another 



1 



PoL 2 



(39) 



to the energy gain due to the tunneling term 

Ej » t\p L 2 , (40) 

where L 2 is the area of each condensate. Note that L 2 
can be significantly larger than the imaging area A. From 
these equations we find that Ej > E c is equivalent to 



t± > 



1 



PoL 2 ' 



(41) 



This condition is compatible with the dominance of A\ 
over A 2 in Eq. ([35]) if 



A> 



1 



PoL 2 



(42) 



In the weak coupling limit A <C p this requirement re- 
duces to A > A/nL A and in the strong coupling 

regime A > fj, Eq. (|42|) reduces to ^jA/nL A < 1. Clearly 
both conditions can be easily satisfied using either small 
imaging areas or systems with sufficiently large number 
of particles nL 2 . 

For convenience we assumed the zero temperature case 
throughout the paper. However, we would like to stress 



that our results will be robust to the effects of tempera- 
ture as long as it is below the temperature correspond- 
ing to the Kosterlitz-Thouless transition. Indeed in cou- 
pled condensates it is the global phase which is destroyed 
by thermal fluctuations, but the relative phase remains 
locked all the way to the Kosterlitz-Thouless transition 45 . 



C. Discussion 

Having the effective action for the fluctuations of the 
order parameter phase we can lift the assumption used 
in derivation of Eq. (f3"2")) that the phases between two 
superfluids are locked together and get 



.4 



A Po A A _ 

.21 w arctan — e 

" 2 u 

x exp 



(e 2 )/2 



Ap A A 
arctan — 



2p t 1 _Ts, D ( A /») L: 



(43) 



We note that ^2,2 is not sensitive to the relative phase 
between the two superfluids and thus Eq. (13"2")) holds for 
an arbitrarily small t±. 
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FIG. 3: Anomalous interference amplitude A2 = ^2,1 + 
A2,2 and separately contributions ^2,1 (Eq. (|43p ) and ^2,2 
(Eq. (132|) ) versus tunneling t±//j, for sample parameters (see 
text for details). The red line (^2,2) is the contribution, which 
comes from normal fermions, the black line is the superfluid 
contribution. 

In Fig. [3] we plot A 2 — A 2 .\ + and separately 
contributions A 2<1 (Eq. (|33j)) and A 2t2 (Eq. p2]0 ver- 
sus tunneling t±/fj, for sample parameters A//i = 0.2, 
A/L 2 = 0.02, TV = 10 3 , where N is the total number of 
fermions per condensate. The graph indicates that there 
is an intermediate tunneling regime, where the superfluid 
contribution A 2 ,i dominates over the normal part ^2,2- 
We note that because ^2,2 is not sensitive to the super- 
fluid gap one can perform separate measurements of A 2 
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in the normal and superfluid regimes. Then the difference 
between the two will be precisely given by A% \ . 

So far considering interlayer coupling we focused only 
on large imaging areas, with both transverse and longi- 
tudinal dimensions large compared to microscopic length 
scales. In the opposite limit one has to study the correc- 
tion to Eq. (U^J) . Instead of repeating rather tedious anal- 
ysis of Sec. IIV Al we will make a couple of simple points, 
(i) There will be additional contribution to 02,2(0) which 
scales as t\. This contribution will have the same ori- 
gin as A 22 in Eq. and will be insensitive to super- 
fluidity, (ii) This unwanted contribution will be greatly 
suppressed because it scales as the square of the imaging 
area, which is small since one length scale is microscopic. 
Thus the effect of the inter-layer coupling on the inter- 
ference in this case will be even smaller than in the case 
of the macroscopic imaging area. 



where we use two imaging areas can not be straightfor- 
wardly applied to bosons because the single particle cor- 
relation functions decay slowly. However, in bosonic sys- 
tems one can deal even without tunneling for the follow- 
ing reason. The relative phase between two independent 
superfluids is random from shot to shot. Nevertheless in 
each shot the interference amplitude is well defined and 
fluctuates only weakly^'- 8 . Thus this unknown relative 
phase can be reliably determined in each run. Then one 
can evaluate anomalous correlation functions putting the 
origin of integration along z in the position of the central 
interference peak. It is easy to see that this procedure 
eliminates the effect of the unknown phase and leads to 
Eq. (USD. 

A. One-dimensional superfluids 



V. ANOMALOUS CORRELATION FUNCTIONS 
IN BOSE SYSTEMS 

So far the main focus of our work was analysis of the 
possibility of measuring anomalous correlation functions 
in fermionic superfluids. There such measurements are 
key for determination of the pairing gap. On the con- 
trary one can get substantial information about the su- 
perfluid properties of Bose systems analyzing normal cor- 
relation functions ^ 18 ' 19 . Nevertheless the possibility to 
measure anomalous correlation functions can provide ad- 
ditional valuable information about properties of these 
systems. As we will see below these functions have 
very unusual behavior in the systems with quasi-long 
order like zero-temperature one-dimensional and finite- 
temperature two-dimensional systems. In this section we 
will give explicit results both for ID and 2D systems. 

We consider setup analogous to that discussed in 
Sec. IIVI Using the same arguments we find that 



A 2 =A d d rF*{v)F 2 {v), 



(44) 



reA 



where F a (r) — (a a (r)a a (0)), a = 1,2; A is the imaging 
area for 2D systems and the imaging length in the ID 
setup. The operators a a have bosonic statistics. As in 
the case with fermions A 2 vanishes if the two systems are 
uncoupled. However, as we argue below if the imaging 
area is smaller than the system size, one can always find 
the regime of small transverse tunnelings t± such that 
the phases of two superfluids are locked together, but 
the correlation functions remain independent of tj_- In 
this case assuming that the two systems are identical, as 
in Sec. ITVl we have 



A 2 ^A d d r\F(r)\ 2 , 

JreA 



(45) 



where F describes anomalous correlation functions in ei- 
ther of the two systems. We note that the scheme I, 



Let us first analyze Eq. (|45|) for two coupled Bose sys- 
tems at zero temperature. If the imaging length is larger 
than the healing length of the condensate then boson- 
boson correlation functions are approximately 



(a(x)a(0))cx(|) 2k ' 



(46) 



where K is the Luttinger parameter related to the inter- 
action strength (see e.g. Ref. [5(|) and L is the length of 
a system. Dots indicate other contributions which scales 
with larger power and thus their contribution is less im- 
portant. For a Lieb-Liniger gas with short range repul- 
sive interactions K 3> 1 in the weakly interacting Gross- 
Pitaevskii regime and K — ► 1 in the fermionized Tonks- 
Girardeau regime. Notice that this correlation function 
increases with the distance. This unusual behavior also 
reflects the scaling of Bq with the imaging length A: 



A 2 (A) ~ CA 2 {Ak)* , Ak < 1 
A 2 (A)~CA 2 , Ak->1, 



(47) 



where C is a nonuniversal constant and k oc \ft±_ is the 
characteristic wavevector corresponding to the transverse 
tunneling (see details below). At sufficiently short dis- 
tances the anomalous interference amplitude A 2 grows 
faster than the area squared. As the imaging size ap- 
proaches the cutoff length the superlinear dependence 
becomes linear and we recover the expected result for 
the system with a long-range order. 

We can make the analysis more quantitative using low- 
energy description of two coupled one-dimensional Bose 
systems (see e.g. Ref. [5l| for more details). In particular, 
for weakly interacting superfluids it can be shown that 
the Lagrangian governing properties of the relative phase 
4>(x, t) between the condensates is: 



47T 



{d x (j>) 2 - 2k 2 cos(0) 



(48) 
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where r is the imaginary time, v s is the sound velocity 
in each condensate, and k 2 rs 4Trb±n/v s K. 

The analysis of either normal or anomalous correla- 
tion functions can be performed using the form-factor 
approach. For more details on this approach we re- 
fer to Ref. for a general scheme and to Ref. [5l| 
for the application to the one-dimensional condensates. 
This analysis is rather involved for generic value of 
the interaction parameter K and requires summation 
of many contributions coming from the intermediate 
processes including creation and annihilation of soliton- 
antisoliton pairs as well as their bound states- breathers. 
Anomalous correlation functions would correspond in 
that case to the soliton-creating or soliton-annihilating 
form- factors 5 ^. On the other hand for weak interactions 
(large K) the Lagrangian can be simplified even further 
if we invoke Gaussian approximation replacing cos(</>) by 
1 — </> 2 /2. More carefully, this can be done using the 
self-consistent harmonic approximation 54 . In this sim- 
ple Gaussian approximation the calculation of anomalous 
correlation function becomes trivial and we get 



F{x) 



a\(x)a 2 (x)a\(0)a 2 (0)) 
2tt 1 

exp 



C(e 



--Ko(kx) 



,(49) 



where C and C are nonuniversal numerical factors and 
A is the short-distance cutoff of the order of inter- 
particle density A ~ n (we assume that k/A <C 1). 
The first exponential factor appearing in the equation 
above is similar to that, which we discussed earlier (see 
Eq. (|4"3"| ). It shows that for small enough t± such that 
K < 2ir/KL the phases between two superfluids are not 
locked and anomalous correlations are exponentially sup- 
pressed. If on the other hand the opposite is true the 
phases are locked together and this term is close to unity. 
And finally the last multiplier gives the spatial depen- 
dence of the correlation function F(x). At kx <C 1 we 
have Kq(kx) r> — ln(Kx) and we recover the asymptotic 
F(x) oc x 2 / K ,in the opposite limit ki » 1 we have 
Ko(kx) <C 1 and thus F(x) « const(a;). 
Using Eq. (gSJ) we find that 



A 2 = AC 



_2_ 


4tt 1 " 


) K exp 


~~KkL 



A 

x J dx exp 
o 



K 



K (kx) 



(50) 



In Fig.|3]we plot the dependence of A 2 /(nL) 2 on kL for 
the imaging size equal to the system size A — L. At small 
tunneling kL <C 1 the two condensates are decoupled and 
the anomalous correlator is exponentially suppressed. As 
the transverse tunneling increases kL ~ 1 the anomalous 
interference amplitude increases faster than linearly and 
at kL > 1 we have A 2 oc L. Note that there is a very 
wide range of parameters 3 < kL < 100, where the su- 
perlinear behavior can be observed. Similarly one can 
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FIG. 4: The ratio of A2 and (kL) vs. kL for imaging length 
equal to the system size A — L; k oc y/t± is the charac- 
teristic wavevector related to the gap A (see Eq. [48}. The 
inset shows A2/A 2 as a function of A/L for different values 
of kL. The superlinear behavior of Aq vs. A is a conse- 
quence of unusual behavior of anomalous correlation function 
F(z)~ (exp[i(cf>(x)+cj>(0))]. 



fix the transverse tunneling and the system size and an- 
alyze the dependence of A 2 as a function of the imaging 
size. The inset shows such plots for different values of kL. 
Again one observes the superlinear behavior of A 2 in a 
wide regime of parameters. We remind that the normal 
interference amplitude Aq defined in Sec. Ill Al always has 
sublinear behavio r 18 i 19 ' 20 . 

We can easily generalize the analysis above to the case 
of finite temperatures. In the regime when the thermal 
coherence length £t ~ Kv s /T is large compared to the 
healing length of the condensates the phase model de- 
scribed the Lagrangian (|4"5| gives the correct low energy 
description of the two coupled superfluids. Within this 
model one finds 



(a(x)a(O)) = nexp 



. ,1 v s \/q 2 +K, 2 

4-7T \ ^ cotn — )L 7yr 9 qx 

2 j — cos^ 



KL 



y/ 'q 2 + ft 2 



(51) 

In the zero temperature limit this expression reduces to 
Eq. (g9]) while at T > v s k we get 



(a(x)a(0)) w nexp 



ttT 
Kkv s 



1 



(52) 



In this case we have exponentially increasing correlations 

'ttT(x - 2/k)" 



(a(x)a(O)} ~ exp 



(53) 



for kx < 1 and then their saturation as kx becomes larger 
than one. In turn this behavior of the correlation func- 
tions implies that in the regime T > kv s the anomalous 
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interference amplitude scales as 



A 2 cx .Aexp — — (A - 2 k) 
\Kv s 



(54) 



at A <C 1/k and then in the usual way A 2 ~ .4 2 in 
the opposite limit. Thus anomalous correlation func- 
tions can be used to probe the temperature in the system. 
We note, however, that in the high temperature regime 
T kv s the anomalous interference amplitude Aq is 
exponentially suppressed. 



B. Two-dimensional superfluids. 

In a similar fashion to the previous section we can 
analyze behavior of the anomalous correlation func- 
tions and the anomalous interference amplitude for two- 
dimensional bosonic superfluids. Below the Kosterlits- 
Thouless phase transition temperature the low energy 
properties of the superfluids can be described by the ef- 
fective Lagrangian (defined as the ratio of the energy den- 
sity to the temperature), which is very similar to (|48|): 



£(r) 



h 2 p s 

AmT 



((V(j)) 2 - 2k 2 cos< 



(55) 



where as before 4> is the local phase difference between the 
two superfluids, T is the temperature, p s is the superfluid 
density, m is the boson's mass, and k 2 w 4mt±/h 2 . Be- 
cause of the formal analogy of Lagrangians (|55|) and (|48|) 
the anomalous correlation function F(r) = (a(r)a(O)) in 
two dimensions is identical to Eq. (|4T))) under the sub- 
stitution K — > ttH 2 p s / \mT) and v s — > 1. So one can 
expect a similar superlinear behavior of the anomalous 
interference amplitude A 2 : 



A 2 ~A 



2+1/2K 



for K\fA *C 1 and A 2 ~ A 2 in the opposite limit. 



(56) 



C. Paired multicomponent bosonic condensates 

Condensates of atom pairs can also be realized with 
bosonic atoms. The original idea of fragmented conden- 
sates goes back to Nozieres and Saint James^. They 
emphasized the difficulty of achieving such states since 
attraction between bosonic atoms favors binding not just 
two but many particles and may lead to the system col- 
lapse. More recently, paired condensates were discussed 
by Kuklov et ali^ for a two component bosonic mixture 
in an optical lattice. Perhaps the most natural setting for 
the appearance of bosonic pairing is polar condensates of 
5=1 atoms in two dimension a 56 ' 57 > 58 . As discussed in 
Refs. 25|26| . in such systems general topological con- 
siderations suggest the appearance of quasi-long range 
order for singlet pairs rather than individual spinor com- 
ponents. This is the system that we focus on in this 
subsection. 



Let tj) a (r) be individual spinor components with m — 
±1,0. We can make a spin singlet pair operator A(r) = 
tp + i(r)ip-i(r) + if} 2 (r). As discussed in Refs. [25U26I ] for 
two dimensional polar condensates, such as 23 Na, one 
expects to find a phase in which (A^ (ri)A(r 2 )) shows 
power law correlations. At the same time there are only 
short range correlations for individual spin components. 
In an interference experiment from a pair of independent 
5 = 1 polar condensates one should measure interference 
amplitude for individual spin components, A m R, then 
construct Ar = A + irA-ir + Aqr. In each shot the 
phase of A(R) is random, so one can again take two re- 
gions, Ri and Rn, and consider (A R Arjj). This expec- 
tation value should decay as a power law of the distance 
between the two regions. 



VI. SUMMARY AND CONCLUSIONS 

In this paper we addressed questions of application 
of interference experiments to detect paired states of ei- 
ther fermions or bosons in low dimensions. We showed 
that direct generalizations of approaches used in ana- 
lyzing interference of independent bosonic condensates 
do not work due to overwhelming shot noise contribu- 
tion. Thus we proposed and analyzed two alternative 
schemes of interference experiments which can be used 
to study anomalous correlation functions, which contain 
information about pairing amplitudes. These (gauge- 
noninvariant) correlation functions provide complimen- 
tary information to normal correlation functions and can 
be used to characterize the properties of the superfluids. 

It was shown how the method of studying the anoma- 
lous functions can be used to detect various pairing or- 
ders. One of the scheme we propose is based on the 
phase sensitive detection employed earlier in the con- 
densed matter systems. On the other hand, another 
scheme deals with two superfluids weakly coupled by in- 
terlayer tunneling. We establish the condition of valid- 
ity of this scheme which involves the tunneling strength, 
imaging area and the system size. We emphasized im- 
portant roles of different form of expansion (transversal 
and longitudinal) and directions of observation. In the 
case of bosonic superfluids anomalous correlation func- 
tions have an unusual property that they increase with 
the separation between quasiparticles. 
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APPENDIX A: EXPECTATION VALUES OF THE 
NORMAL CORRELATIONS FOR THE SYSTEMS 
WITH PAIRING 

In the BCS approximation the expectation value of the 
square of the interference fringe amplitude discussed in 
Sec. II can be found as 



the interactions are not important during the expansion 
then the spatial image of the cloud after time of flight 
gives the momentum distribution of atoms in the initial 
condensate. As it is shown in Ref. @ the density-density 
correlation maps to the Fourier transform of the pairing 
amplitude: 



(\A\ 2 )=A 



(2tt)° 



in 



(Al) 



(n(R,i)n(-R,t)) = |f(Q) 2 



(Bl) 



where nf(Ejs) is the Fermi distribution function, E^ = 
Y^ + AJ, £k = k 2 /2m-/i, A is the imaging area of the 
interference. The pairing gap function A& is a constant 
Aq for the s-wave pairing and the fc-dependent function 

Afe = A^\cos(k x ) — cos(fc !/ )) for the case of the d-wave 
pairing. In the absence of superconductivity and at zero 
temperature \A\ 2 = 2N, where N is the total number of 
fermions in each system. The nonzero value of |j4| 2 in 
this case simply reflects anti-bunching of fermions (note 
the negative sign in Eq. ((3j). As expected for fermions 
A~y/N (see e.g. Ref. [if). 

Considering s-wave pairing and using the constant den- 
sity of states po in 2D one finds 



| A | 2 = 2N-AA ( Q s) p (~ + arctan(i s ) - 2[t 2 



where t s = ///Aq and p is the chemical potential. In the 
BCS limit p 3> A^ the equation above reduces to A 2 ~ 
2N — AirpoA^ . If we formally extrapolate BCS theory 

towards the unitarity limit p — + and Aq S ^ — ► n/po then 
A2 — > N(2 — 7r/2). So we see that A2 is a monotonically 
decreasing function of the pairing gap. The results will be 
somewhat different for the superfluids with the d— wave 
pairing in which case is 



(A2) 



(\A\ 2 ) = 2N + 2Ap A^\[t 2 d + l 



(A3) 



2tE(-^ 2 ) + 2(l + i 2 )2 E ( 



TT7f» ^ 



where t d — p/ A^\ 



APPENDIX B: TRANSVERSAL VS. 
LONGITUDINAL EXPANSIONS 

In this Appendix we compare the two regimes of expan- 
sions: the transversal expansion advocated in the main 
text versus longitudinal one which is more close to to the 
standard time-of-flight technique. 

The quantity similar to A2 in Eq. (|12[) can be also mea- 
sured in the standard time of flights experiments. For the 
low-dimensional superfluidity it is advantageous to have 
longitudinal expansion so that the atoms from different 
layers do not mix with each other. If one assumes that 



where Q = mR/ht. We note that the transverse ex- 
pansion discussed in the previous section directly probes 
the spatial structure of (the square of) the pairing wave 
function and thus gives a complimentary information to 
the quantity (|B1|) . If one integrates Eq. (jBip over the 
momentum Q then one recovers the expression identical 
to Eq. (TIB")) . So the two setups give equivalent informa- 
tion. We note, however, that in the transverse expansion 
regime one can benefit from independently averaging over 
several imaging areas within a single shot. 

Doing the longitudinal expansion one can also deter- 
mine the spatial structure of the pairing amplitude. In- 
deed if one integrates Eq. (|B1[) along a preferred direction 
then the outcome should be isotropic for the s— wave 
pairing and highly unisotropic in the d— wave regime. 
Thus within the BCS model in the s— wave case and in 
the d— wave case if one integrates along the anti-nodal 
direction (where the pairing gap is maximal) one gets 



J dx(n(R,t)n(-K,t)) oc A 



p+ y/pJ+A 2 

vV + A 2 



(B2) 



On the other hand if one integrates Eq. (|Bip along the 
nodal direction in the d— wave case (i.e. the direction 
where the pairing gap vanishes) one should get zero. 

The clear advantage of the longitudinal expansion 
method is that one avoids the issues with coupling be- 
tween different layers. In principle, one can perform this 
experiment even on a single layer. However, there is a 
big disadvantage too. Namely, one has to rely on the 
free expansion to get the right correlation functions. In 
the case of a bilayer system predominantly transverse 
expansion is guaranteed by the large kinetic energy of 
the transverse confinement. On the other hand for the 
longitudinal expansion collisions between atoms shortly 
after expansion started can affect the outcome of the ex- 
periment. We note that instead of purely longitudinal 
expansion one can have the full three-dimensional time 
of flight experiment and get similar results. However, if 
one is interested in two-dimensional superfluidity allow- 
ing expansion in all three dimensions one clearly looses 
in the signal to noise ratio. 
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APPENDIX C: EFFECTIVE ACTION 
APPROACH TO COUPLED SYSTEMS 

The effective action (S e g = J °° £(t)c?t) for the phase 
fluctuations reads^: 



n=l 

~ J d<W(|Ai(r,r)| 2 + |A 2 (r,r)| 2 ), (CI) 

where g is the interaction strength, which we assume to 
be short range, Go and £ are 4x4 matrices: 



So 



d 
G, 



Si T 
Tt E 2 



(C2) 



Here Gi,2 are the fermion's Green's functions in the su- 
perfluid with a non-fluctuating phase. In the Nambu 
notation their inverses are: 



^1,2 



|A ll2 | a T + A_v 2 



(C3) 



Ei t 2 contains fluctuations of the order parameter—: 

h 2 

Ei >a = -^(V 2 6»!, 2 + 2V6»i, 2 V)t 
4m 

^1,2 + |-(V0i, 2 ) 2> ) r z - 8\A lt2 \T x . (C4) 
2 8m / 

And finally T corresponds to the tunneling coupling be- 
tween two layers: 



T = t ±c i(ei/2 ~ e2/2)T ^ 



(C5) 



For simplicity we assume that two layers are identical 
and | Ai| = |A 2 | = A 

Next we expand the effective action (|C1[) to the second 
order in small fluctuations in derivatives of $i j2 and to the 
second order in t±. For simplicity we ignore fluctuations 
in the magnitude of the order parameter <5|Ai j2 |. Then 
one obtains: 

S cS « Si + S 2 + S 1<2 , (C6) 

where Si j2 correspond to the quadratic Lagrangians 
of decoupled layers and 



5ia = -Tr(G TG T t ). 



(C7) 



Note that we can ignore slow spatial variations of the 
phases #i. 2 in the tunneling matrix T in the equation 
above. Indeed keeping these variations will result to cor- 
rections to the gradient terms in Si ;2 proportional to 
tj_ • Then taking into account the expression for the non- 
diagonal elements of the Green function, 



Go — — 



iuJT + g fc r 3 - A fc Ti 



(C8) 



where £fc is dispersion and 10,1,3 are Pauli matrices in 
Nambu space, Eq. (|C7|) can be evaluated in the limit of 
zero temperature as 

Si 2 = — %t±Po 

dw f d 2 k 
2^ 



X C0S(6>i 



r-(C9) 



(2tt) 2 (lo 2 + El) 2 ■ 

where uf(0) is a density of state at the Fermi level and 
= + A 2 .. For the s-wave superfluid we use A& = 
A and = k 2 /2m — /1 whereas for the d-w&ve pairing 
A fe = Ao(cos(fc x ) - cos(fcy)). 



APPENDIX D: TECHNICAL DETAILS 

In this Appendix we present some details of our compu- 
tations of contributions from the anomalous correlation 
functions. While the expression for the most functions in 
the case of s-wave pairing can be evaluated straightfor- 
wardly, its treatment for the d-wave case require some ap- 
proximations which we describe here. However, some in- 
tegrals can be done without these approximations which 
is also discussed. 



1. Nodal approximation 

In the T — » limit we are interested in, to evaluate the 
d-wave related integrals we use the nodal approximation 
developed in£2. It consists in focusing on the regions close 
to four nodes of the d-wave order parameter, k„ = kpk n , 
(n = 1, 2, 3, 4) on the Fermi surface. Explicitly, 



ki 



V2 ' 
-x-y 

V2 ' 



k 4 = 



V2 ' 
x-y 



V2 



(Dl) 
(D2) 



In the vicinity of nodes we can expand k = k„ + (5k, and 
then 



£fc = v F 8k±, A k = w A <Sfc||, (D3) 
v 2 F Sk 2 ± + v 2 A 5k 2 (D4) 



El 



Here <5fc_i_,|| ar e components of momentum perpendicular 
and parallel to the Fermi surface, up is a Fermi velocity 
and va = |9Ak/9k|||. Finally introducing the angular 
parametrization 

£k = vpSk± = ecosa, A& = VA.8hn = esina(D5) 
E k = e (D6) 

The integrals can be performed according to the following 
rule: 



ri 2 k 

(27T) 2 
1 

2ttvfVa 



(...) 



4 

£ 

n=l 



ddk\ddk 



(2tt) 2 



(...) 



4 



(D7) 
(D8) 
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where the limit e max ~ Ao makes sure that the integra- 
tion area is equal to the area of the Brillouin zone. It can 
be send to infinity for practical purposes. 

Considering the tunneling integral of Eq. (|C9[) we ob- 
serve that 



2tt 2 



VfVa 



du) 
2^ 



arctan(y) — y log 



kdk 



1 



(D9) 
(D10) 



where y = uJmax/tmax- In the limit y — > oo the expres- 
sion in square brackets goes to tt and we arrive to the 
equation (|38|) if one assumes that the ratio e m ax/vFi>A is 
a constant independent of A. 



2. Explicit evaluation of some integrals 




FIG. 5: Comparison of the tunneling elements for the s-wave 
pairing (lower curve) vs d-wave pairing (upper curve) as a 
function of 

Some explicit expressions are possible to obtain for sev- 
eral important quantities: in particular we focus on the 
function 6q introduced in Sec. IIII CI and the tunneling 
integrals. 

Starting from the expression for the anomalous Green 
function 



for x = A 2 cos(40)/2 one can obtain a systematic ex- 
pansion of the integrand. Restriction to the first term 
in this expression already produces a very good uniform 
approximation for the integral. Doing the <f> integration 



jkf cos(0) cos ( 2 0) = -2tt J 2 (kf) 



where f = r(2m* /h 2 )z gives an expression 



F(r) w -p A / kdk 
Jo 



J 2 (kr) 



sj(k 2 - + 



A2 

2 



(D13) 



(D14) 



where po = m* /irh 2 is the constant density of states 
in 2D (summed over spin polarizations). The expres- 
sion (|D14p can be further evaluated using the so-called 
^-approximation, frequently used in the theory of super- 
conductivity. In the argument of the Bessel function we 
write k 2 — p = (k — yfp)(k + ^fp) w (fc — ^fp)2^fp. On 
the other hand, in the BCS regime where p/A >> 1 
after variable's rescaling we shift the lower limit of inte- 
gration to — oo and use the summation theorem for the 
Bessel function J v (z + t) = Y^n=-oo Jv~k{t)Jk{z). Re- 
stricting to the n = from this sum already produce a 
good approximation for the function F(r), 

W=-^MVMlo^)Ko{^) (D15) 

Now, the integral for 6g can be computed in the same 
BCS limit by using the asymptotic behavior of Bessel 
functions for large and small arguments. The log 3 -result 
is similar to the s-wave case but with different prefactor 
given in Eq. ([22]) . 

Tunneling integrals in the rf-wave case can be evalu- 
ated exactly. The answer is shown in Eq. ([38]) with the 
function 



2 , n y 2 



T D (y) = 1 + ^[ V / TT^E(-^) 



iry- 1 + y 2 

E(-y 2 ) K(-y 2 ) 



(D16) 



2^ vl J K(-f_)] 



^T+y 2 i + r 



F(r) = 



2m* kdk 



2ir ^ikr cos < 



A cos 2(j)d(j) 



o (2Tr) 2 h 2 J y/(k 2 -p) 2 + A 2 cos 2 20 

(Dll) 



writing cos 2 20 = (1/2) + cos(40>)/2 and using the expan- 
sion formula 



v 7 ! + x 



{~l) k {2k)\x k 

k=0 



(k\) 2 2 2k 



(D12) 



where E(y) and K(y) are complete elliptic integrals and 
y = A/ p. In Fig. (|5|) it is compared with exact expression 
for the s-wave tunneling integral 



Ts(y) = 1 



l 



(D17) 



The similarity of the tunneling integrals for the s— and 
d— wave pairings is remarkable. 
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